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Abstract 

The densities of states are evaluated for very short chain molecules made up of overlapping 
monomers, using a model which has previously been shown to produce helical structure. The results 
of numerical calculations are presented for tetramers and pentamers. We show that these models 
demonstrate behaviors relevant to the behaviors seen in longer, helix forming chains, particularly, 
"magic numbers" of the overlap parameter where the derivatives of the densities of states change 
discontinuously, and a region of bimodal energy probability distributions, reminiscent of a first 
order phase transition in a bulk system. 
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I. INTRODUCTION 



Helices are a common structural motif in biological molecules, from the a-helix in proteins 
to the culturally iconic double helix observed in double-stranded DNA. In a living cell, the 
adoption of stable helical structures allows these molecules to place functional groups in 
specific positions and orientations, and holds the polymer backbone away from the solvent, 
protecting it from chemical attack. The consensus view of helix formation follows the work 
of Pauling et al. jl|; biological helices are stabilized by orientationally-dependent hydrogen 
bonding, with their chirality arising from the chirality of the polymer molecule. 

Those same properties which make helical molecules so useful in living cells also make 
them useful in the context of nanotechnology. Unfortunately, while our understanding of 
biological helices is good at explaining why polypeptides and polynucleotides do form helices, 
it does not provide useful prescriptions for developing alternative helix-forming molecular 
architectures. To gain the understanding necessarv to develop such prescriptions, many 
workers have considered "reduced models" for helix formation, which attempt 

to capture the underlying physics of the phenomenon in as simple a manner as possible. 

Over recent years, simulation studies of such reduced models have yielded surprising 
results. In particular, several polymer models have been proposed which produce helical 
structure while interacting via isotropic potentials ; that is, helix formation without 

"designed-in" preferred interactions, with spontaneous chiral symmetry breaking. Maritan 
et al. have shown that helices are "maximally compact" structures for string-like 

objects. This suggests that helix formation arises from geometric symmetry breaking, akin 
to crystallization. A better understanding of how this symmetry breaking can arise should 
lead to the better prescriptions for helix-forming architectures. 

In the study of a helix formation in polypeptides, the starting point is the observation 
that helices are quasi one dimensional objects, which can be looked at as a spin chain. The 
standard approach jjj], 10, H] is to attribute amino acid residue conformations to spins, 



either H ("helix", that is, capable of forming a hydrogen bond compatible with a helical 
structure) or C ("coil" , otherwise). A spin chain representation is then made up of these 
states; in the simplest form [3], residues which are neighbors along the peptide backbone 
interact according only to their H/C attribution and amino acid type. Modern versions of 
this approach Q| include many-body "capping interactions", which are non-pairwise, non- 



2 



local interactions between residues; the strength of these interactions, however, still depends 
only upon the residue type and H/C attribution. Such models have achieved considerable 
success in helical structure prediction for polypeptides. For more general helix-forming 
systems, the proper attribution of a backbone segment to "H" or "C" type is not clear. 
However, the success of the spin chain approach to helix formation in polypeptides suggest 
that a similar approach may be fruitful. 

For a linear polymer of spherically symmetric monomers, single monomers are not the 
equivalent of amino acid residues for helix formation, as they have no internal degrees of 
freedom. From symmetry arguments, the minimum possible such building block must be a 
tetramer; helices break chiral symmetry, and a tetramer is the shortest length chain which 
may exhibit chirality. Similarly, the behavior of a pentamer should contain information on 
how neighboring chiral centers interact, and so forth for longer chains. 

In this paper, we seek complete enumeration of the partition function for tetramers and 
pentamers, using a simple polymer model which has previously been shown to produce 
helices 3L This enumeration is performed using a methodology similar to that followed by 
Taylor [12|| for short tangent square-well chains. The intention is to identify the building 
blocks necessary for helix formation in longer chains, and the origins of the behaviors which 
allow helix formation in longer chains. The methodology and results of this enumeration 
are intended as a staging post for the construction of generic spin-chain models of helix 
formation 

The remainder of the paper is structured as follows. In Sec. [Ill the polymer model which 
is to be studied is described. In Sec. IHIl the method by which the partition functions for the 
model are calculated is described. The results calculated from these partition functions are 
described in Sec. HVl Finally, in Sec. [V], these results and their implications are discussed. 



II. MODEL 



The polymer model consists of a linear chain, bond length /, of N hard spherical monomers 
with diameter a. The degree of overlap between monomers is determined by the reduced 
parameter a/l. For a/l = 1, this is the familiar tangent sphere polymer model. We consider 
chains with a/l > 1, that is, with overlapping monomers. Interactions between non-bonded 
monomers (separation r) are given by an isotropic square-well potential: 
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Figure 1: Schematic phase diagram from simulation for a helix-forming 20mer, as described in the 

f. fl. 



main text. Reproduced from Ref. 



oo r < a 

u(r) = { -e a <r < Xa (1) 
Xa <r 

where A is the well width (taken as 1.5 in this work), and the well depth e sets the energy 
(and hence temperature) scale. We follow the protein literature, by denoting interactions 
between particles where a < r < Xa as contacts, and interactions where r < a as overlaps. 
Interactions between monomers separated by two bonds along the chain are referred to as 
1-3 interactions; interactions for monomers separated by three bonds are referred to as 1-4 
interactions, and so forth. 

In previous simulation work, we have used a version of this model where individual bond 
lengths were allowed to vary by ±10%. It has been suggested that such bond length varia- 



tion can enhance the ergodicity of a simulation compared to rigid bonds [13| ; further, this 
allows the configurational and momentum parts of the partition function to be factorized. 
With such bond length fluctuation, the system has been shown to form helices for 20mers 
(polymers of length N = 20). The observed phase diagram is shown in Fig. [U the system 
is observed to form two distinct helical phase, "helix 1" (stable at higher temperatures, and 
with a smaller radius) and "helix 2" (stable at lower temperatures, and with a larger radius). 
The following work does not include such bond flexibility, as the extra degree of freedom 
per bond would make the problem very much less tractable. 
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Figure 2: A core-softened potential. If the shoulder height A S> kT, the effective core diameter will 
be <7 rather than I; hence a chain of monomers interacting via such a potential with bond length I 
would act as the overlapping square well monomer model presented here. 

A. Physical Relevance 

With any such "reduced model", however interesting the behaviors, the question of phys- 
ical relevance must be answered. The idea of overlapping monomers is consistent with the 



Van der Waals radii of atoms in "realistic" potentials such as CHARMM [14|], where atomic 
radii are often larger than the bond length to neighboring atoms. On a larger scale of ap- 
proximation, if amino acid residues are approximated by interacting spheres, the radii of 
gyration for amino acids can be larger than their center of mass spacing along the peptide 
chain. 

Since protein molecules form intra-chain hydrogen bonds, an interesting parallel can be 
made to "core-softened potentials" (see Fig.EJ), which have been used to study the anomalous 



behavior of water 
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12|| . These isotropic potentials have a shoulder (diameter a) around 



a repulsive core (diameter /), representing close packed but non-hydrogen bonded pairs, and 
an outer well (diameter Act) which represents hydrogen bonding interactions. In a chain 
of such monomers with bond length /, if the difference between the potential energy of 
the shoulder and the potential energy in the minimum is sufficiently larger than &bT, the 
effective repulsive core diameter will be the shoulder diameter; at low temperatures, a chain 
of such monomers would behave as the overlapping square well monomer model presented 
here, forming helical structure. 
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Figure 3: A cartoon of the model. Five monomers are shown, of diameter a, bond length I, with 
bond angle 6\ and the two dihedral angles <f>\ and indicated. 

III. METHODS 

We consider 4- and 5-length polymers of the type described above, as shown in Fig. [3l 
The position of monomer i is denoted by R,. Bond vectors are defined as r-j = Rj + i — Rj. 
Separation between monomers i and j is denoted r^ . The bond angle around monomer i 
is defined as the angle between bond r^_i and r^, that is cos#j = • i-j. The dihedral 
(torsional) angle <pi is defined as the angle between the planes formed by the vector pairs 
(r^-i, Ti) and (r*, relative to the cis conformation (i.e. <fi = is cis, <fi = n is trans). 

We take positive as a right-handed rotation. This is, however, arbitrary, as the underlying 
model is achiral. We do not consider translations and rotations of the entire molecule; as 
such, we fix the position of the first monomer, as well as the plane made by the vectors r x 
and r 2 . The configurational integral of such an n-mer is given by Z n , defined as: 



'n— 3 „,r \ in— 2 



Zn = (il/ III/ f^QideA x 



J=l"-« / \j=l 

exp(-/3£({R})) (2) 

where E ({R«}) is the total configurational energy for the system (the sum of the pairwise 
interactions as Eq. CD) and j3 = 1/ksT, the inverse temperature. 



A. Tetramer 



We initially consider a tetramer. The configurational integral is given by: 
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Z A = I I d(f>i / sin^irf^i / sin 9 2 d9 2 x 

J-TT Jo Jo 

exp (-(3 (u (r 13 ) + u (r 24 ) + u (r M ))) (3) 



Separations are given by: 



r\_ w = 2/ 2 (l + cos^) 

/ Xi 



(4) 



and 



\2 



r i-i,i+2 (0i,8i+i,<i>i) = I 2 ((I + cos 9i + cos 9 l+1 y 

+ sin 2 9i + sin 2 9 i+1 
—2 sin 9i sin 9 i+1 cos 4>i) 
= l 2 ((x i + x i+1 -2) 2 /4 

(4 - ^) /4 (5) 
+x i+1 (4 - /4 



- cos <pi^/xi (4 - Xijx 

Since we are working with variables of squared separation, for notational convenience we 
also define a = (cr/l) 2 . Physical bounds for Xi are a < Xi < 4, since a separation of less than 
a represents an overlap. It is natural to switch variables in the configurational integral to 
separations x,, giving: 



c4 

*C2 X 



/7T /*4 /*4 

#! / darx / dx 2 
7T JO JO 

exp — (3(u (ly/Xi) + U (/v^) 

+u (ly/yi (xi,x 2 ,0i)))) (6) 

Since we are working in a square- well system with discretized energies, it is now convenient 
to switch to a density-of-states representation: 

3 

Z 4 = (l/4)$>4(£0exp(/3eA;) (7) 

k=0 



where uj n {k) is the density of states for the n-mer with k contacts. For the tetramer, we can 
write the appropriate integrals: 



d(f)i / dxi / dx 2 x 

J min(4,A 2 a) J min(4,A 2 a) 

(y x (xi,x 2 , 0i) - A 2 a) (8) 



/7r />4 /•min(4,A 2 a) 

C?01 / dX\ \ dX2 X 

-7r •/ min(4,A 2 a) »/ a 

(yi (xi,x 2 ,0i) - A 2 a) 

/it ;>4 />4 

c?0i / dx\ \ dx 2 X 

-7T J min(4,A 2 a) J min(4,A 2 a) 

[0 (yi (x 1 ,x 2 ,(f>i) - a) 



-0 (yi(x u x 2 ,(j) 1 ) - \ 2 a)] (9) 



/7r pi ;>min(4,A 2 a) 

dcpi I dxi I dx 2 X 

■it »/ min(4,A 2 a) J a 

[0 {xi,x 2 ,<p 1 ) - a) 

- (y 1 (x 1 ,x 2 ,(p 1 ) - \ 2 a)] 

pi: />min(4,A 2 a) /•min(4,A 2 a) 

+ / c?0i / dx\ \ dx 2 x 



7r */ a 



©(^(x!,^,^) - A 2 a) (10) 



/7r /•min(4,A 2 a) /•min(4,A 2 a) 

c?0i / <ixi / dx 2 x 

it J a J a 

[0 (yi (xi,x 2 ,0i) - a) 

-e(y(x 1 ,x 2 ,<l> 1 )-X 2 a)] (11) 

where (x) is the Heaviside step function. The min terms exist to prevent unphysical limits 
of integration when A 2 a > 4 (in which case 1-3 contacts are "always on"). We note that 
uiiik) is a sum of integrals of the general form: 

U = 2 d(j)i I dxi J dx 2 Q(yi(x 1 ,x 2 (f) 1 ) - h) (12) 

JO Jxu Jx 2 , 



where we have used the symmetry of the system to simplify the integral. The integrand 
is non-zero for that region of (xi,x 2 ,(pi) space for which y\ > h. We can solve Eq. to 
find the bound of this space with respect to Xi (or, by symmetry, x\) for given ((f), h), which 
we call x c : 



x c (x, (j), h) = 2 yx (t (4-x) + (h- 1)) 

+sgn(cos(0i)) a/ t (4 — x) x x 
yj (2x (2t + h-l)- tx 2 -(h- l) 2 )^) 
I [x (x + t(4-x))) 



(13) 



where we use t = cos 2 0, and sgn(x) returns the sign of x. Similarly, we also solve for the 
value of at which = h for given X{ and Xj+i, which we denote </> c : 



4> c (xi,x i+ i,h) = < 



7T $ < -1 

arccos($) |<&| < 1 
$ > 1 



(14) 



where the ratio <E> is defined as: 



2 + XiX i+ i — 2h 



y/x~i (4 - Xi) x i+1 (4 - x i+1 ) 
We consider the shape of this boundary: 



(15) 



dxi 



Xi + 1 



2a/^+i (4 - x i+1 ) 
sin (pc^JXi (4 - x^ 

(2(1-/0-2* ( Xm + (l-/i))) 



(16) 



The only part of this equation which can be negative is the final bracket. Within the 



range a < £j < 4, a < h < A 2 a, 1 < a < 4 and A 2 > 1, it can be shown that 



is 



dx c 



non-positive. By symmetry, Q 
([H as: 



is therefore also non-positive. As such, we can write Eq. 



U = 2 



i" rxih rx 2 h 

d(pi / dx\ 



dx? 



min(z 2 h,max(x 2 ;,a; c (a:i,(/i,/i))) 



(17) 
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The max term picks the larger of the original lower limit (X21) and the value of %2 below 
which the Heaviside function integrand in Eq. (TT2l) becomes zero. The min prevents the 
unphysical result of the lower limit becoming larger than than the upper limit. 

We note that, by symmetry, if x\ = x c (x2,(f),h), then X2 = x c (xi,(f>,h). Hence, we 
can immediately see that the solution to X2h = x c (xi,fa,h) with respect to x\ is X\ = 
x c (%2h> 0! h). For x\ < x c (x2h, 0, h), the upper and lower limits on the innermost integral 
are equal, and hence the contribution to the integral is zero, hence: 



f 4 — 2 I dfa I dxi x 

J min(x lh ,max(x 1 i,x c (x 2 h,<t>,h))) 

rX2h 

dx 2 (18) 

' max(x 2 i,x c (x 1 ,(f>,h)) 

Since x c (x, 0, h) is a monotonically decreasing function of x across the range of interest, 
we can now propagate the max term in the middle integral out to the dihedral integral: 



/•t r x ih 
/ 4 = 2 / dfa / dx\ x 

J 4> c {x u ,x 2h ,h) Jx u 

/ dx 2 

J max(x2i,x c (xi,(j>,h)) 

r<Pc(x ll ,x 2 h,h) rxih 
+2 / dfa I dx\ x 

JO Jmin(x lh ,x c (x 2h ,(f>,h)) 
x 2 h 

dx 2 (19) 

max(x 2 i,x c (xi,(f),h)) 

The remaining max and min terms can then be propagated out in a similar fashion, remem- 
bering that C (x, x', h) is a monotonically decreasing function of x and x' across the range 
of interest: 



/ 4 = 2[;, , dfa r xh dxx r 2h dx 2 

J* J4>c(xu,x 2l ) ^ L Jx u L Jx 2l z 

+2 f^2? l) ^ < « dfa fj , dxt f X2h dx 2 

Jma,x(<j> c {xii,x 2h ),(t) c {x lh ,x 2 i)) ^ 1 ->x c {x 2l ) 1 J x 2l 

+2 f^(*»f») dfa f Xc{x2l) dx! f X2h dx 2 

Jmax(4> c (xu,x 2h ),4> c (x lh ,x 2l )) t i J Xll 1 Jx c £ 

r^ c (xu,x 2h U c {x lh ,x 2l) ) rx lh rx 2h ^0) 

J<f>c(x U ,X 2h ) ^ LJ X U 1 JX C Z V ' 

+2 f^f" , . , . dfa f x ) h , dx x f X2h dx 2 

Jmm(<p c (x ll ,x 2h ),4> c (x lh ,x 2 i)) ^ 1 •> x c (x 2 i) 1 Jx 2t 

2 fM^2 h ) d<j) r*M d p h d 

Jmm(4> c (x u ,x 2h ),4> c (x lh ,X2i)) ^ L J x c {x 2h ) 1 Jx c z 

™(Wx 1| Ia )A( Ilk a)) d(j) d r*2H d 

J<pc(x lh ,X 2h ) ^ 1 Jx c {X 2 h) «c 
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Figure 4: Graphical interpretation of the integrals in Eq. (|20l) when <p c (xih, x 2 i) > 4> c (xihX2h)- 
Part (a) shows the entire space in x%, x%, <f>i; the total integral is the volume above the gray surface. 
Part (b) breaks the volume corresponding to the total integral into the parts listed in Eq. ([201) 

For notational convenience, superfluous arguments to the functions <p c and x c have been 
omitted; that is, h and variables of integration. There are min and max terms in the dihedral 
integral since, without knowing more about the original limits Xu and it is not possible 
to tell whether C (xi h , x 2 i, h) < <f) c {xu,x 2 h, I)- 

The first integral is trivial. The remaining integrals include a term J x c (x\i, 0i, h)d<f>i) 
this is analytically tractable, resulting in terms involving elliptic integrals, but is simpler to 
treat numerically. The final two terms include integrals of the form J J x c (xi, <pi, h)dxid(j), 
which are not analytically tractable, and are thus treated numerically. 

Calculation of the f± integrals via Eq. f l20l) allows calculation of the density of states 
via Eqs. ([8] -ELD), from which the equation of state and energy probability distributions 
can be determined. Structural information, in the form of the dihedral angle probability 
distributions, is also easily accessible. The dihedral density of states g^i, k) is given by 
integrals as Eqs. j8] - EEE]) without the integral over the dihedral angle. For the tetramer, 
these give tractable though lengthy analytic forms. The probability of observing a given 
dihedral angle is then given by: 



4 




(21) 



fc=0 
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B. Pentamer 



An equivalent procedure may be carried out for a pentamer. We first introduce the 1-5 
separation, z t (xi, x i+1 ,x i+2 , (f>i, 4> i+ i) = r?_ li+3 /l 2 : 

Zi = Xi + x i+2 - yxi+i ( 4 - x i+l) x 

cos(j) i+1 Xi^x i+ 2 (4 - Xi +2 f) /4 

+ (X i+ i - 2) ^XiX i+ 2 - COS (f>i COS (f> i+1 X 

\Jxi (4 - Xj) x i+2 (4 - x i+2 )) /4 
+ sin 0j sin x 

yjxi (4-x l )x i+2 (A-x t+2 )/2 (22) 

where the arguments to Zi have been omitted. The partition function for the pentamer is 
given by: 

6 

Z 5 = (1/8) ^{k) exp{(3ek) (23) 

Equivalent expressions to Eqs. (IMTll are simple to construct, using the equivalent form 
to Eq. (fT2l) : 




e(y 1 -h 1 )e(y 2 -h 2 )e(z 1 -h 3 ) (24) 

where we have suppressed the arguments of yi, y2 and z\ for notational ease. This integral 
is constructed (without loss of generality) such that 4>i i s always right-handed. Explicit 
bounds of integration due to 1-4 interactions can be treated in the same manner as for the 
tetramer case. Bounds for the X\ and x 2 integrals as a function of 4>i, and for the x 2 and x 3 
integrals as a function of 02, are determined exactly as Eq. (|20| ). This leads to single ranges 
of integration for X\ and 23, and two sets of ranges of integration for x 2 . The proper range of 
integration over 22 is then the overlap of these two ranges. Explicitly treating the bounds of 
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Table I: Comparison of the densities of states for a square well tetramer chain, a/l = 1.0 and 
Act = 1.5 calculated in this work (uj^(k)/Al , with the factor 1/4 (as described in section fTTT|) and 
by Taylor [3] (uj^^^^k) , or in the original terminology). Suppression of the unimportant 
multiplicative factor of 8ir in the work of Taylor leads to the difference in the values; it can be seen 
that including this factor, the values differ only in the fourth and fifth significant figure. 

k u 4 (k)/4l 4 uf aylor \k) w 4 (fc)/ (327rJ^ aylor \k) 



4.78131 0.19029 

1 5.59121 0.22247 

2 2.42528 0.09650 

3 0.62013 0.02467 



0.999750 
0.999987 
0.999986 
1.000170 



integration due to 1-5 interactions is not trivial, and as such the resulting integral is treated 
numerically. Dihedral probability distributions P (0i, 02, T) can be calculated from dihedral 
densities of states #5 (0i, 02, k) in an analogous manner to the tetramer. 

IV. RESULTS 

Using the results presented in Sec. IHIl we have evaluated the full partition functions for 
tetramers and pentamers. Results for the tetramer have been calculated with the Mathe- 
matica symbolic algebra package, using Gauss-Kronrod numerical integration. Results for 
the pentamer have been calculated using ten-point Gauss-Legendre quadrature Both 
methods of integration have been checked by comparison against the tangent chain results 



presented by Taylor 



12|. The pentamer results have been verified against short Monte Carlo 



simulations (data not shown). 



Tetramer 



To validate the method, we compare our calculated densities of states for tetramer tangent 
square well chains (a/l = 1) to those presented by Taylor [12||. These results are shown in 
Table CO It can be seen that the results are equivalent to four significant figures aside from 
an unimportant multiplicative factor. The method of Taylor does not use explicit limits of 
integration, instead numerically integrating the Heaviside functions in Eq. (fT2lk strictly, the 
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Figure 5: Densities of states for tetramers plotted against a/l. 

method presented here should be more accurate, though these results suggest the difference 
is not significant. 

The calculated densities of states as a function of a/l are shown in Fig. El The densities 
of states for the k = and k = 1 states are zero for a/l > 4/3. For overlaps greater than 
this "magic number", 1-3 interactions become "always on" — that is, fj-i^+i < Act for all 
values of Q{ with A = 3/2. This also gives rise to a kink (discontinuity in the derivative) 
of 074(2,0"//). This is because at a/l = 4/3, the first term in Eq. ffTOl) (which refers to the 
density of states for tetramers with a single 1-3 contact and a 1-4 contact) becomes zero, as 
the limits on the x\ integral become equal. 

As a/l — > 2, the polymer becomes increasingly rigid, and the available conformational 
space vanishes. The calculated densities of states show the correct behavior at this limit. 

Properties calculated from these densities of states are shown in Figs. [6] (energy) and [7] 
(heat capacity). 

At a/l = 4/3, the slope of the energy and heat capacity contours show discontinuities in 
their derivatives. As such, the derivatives (dU/da) T and (dCy / 'da) T have singularities at 
a/l = 4/3, however, these are not physically meaningful response functions. In simulated 
systems [3j|, bond lengths are not rigid, and bond length fluctuations will have the effect of 
"smoothing out" the discontinuity. 

Though the tetramer does not show any other discontinuities, it does show a line of 



maxima in heat capacity with respect to temperature. We follow Taylor 



121 ] and Zhou, et 



al. [13j| in ascribing these maxima to collapse of the tetramer into compact conformations. 
The strength of these maxima can be seen to decrease with increase in a/l. Further, the 
line of maxima shows re-entrance with respect to a/l, with the temperature at which heat 
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Figure 6: Ensemble average energies (E) /e for tetramers plotted against a/l and temperature T. 
Solid lines show energy contours at the labeled value. Note the discontinuities in the slope of the 
energy contours at a/l = 4/3 (shown by the dotted line). 
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Figure 7: Configurational heat capacity Cy for tetramers plotted against a/l and temperature. 
Solid lines show contours at the labeled value. The dashed line shows the line of maxima in Cy 
with respect to temperature T . Note the discontinuities in the slope of the heat capacity contours 
at a/l = 4/3 (shown by the dotted line). 

capacity is maximal itself having a maximum with respect to a/l at a point below a/l = 4/3. 

Representative results for the torsional behavior of the tetramer are shown in Fig. El 
where we show the probability P(<pi] T) for four values of a/l. At low values of the overlap 
(a/l < 1-48), we see maxima in P{4>\\ T) for non-zero 0i at all temperatures, with the 
maxima becoming stronger and moving closer to zero (cis conformation) as temperature 
decreases. For intermediate values of overlap (1.48 < a/l < y^(3 + v / 5)/2), weak maxima 
in P(0i;T) are seen for non-zero (f>i at high temperatures; however, the most probable 
conformation becomes (pi = (cis conformation) at low temperature. For large values 
of overlap (a/l > J (3 + v / 5)/2), P(4>i\T) has only a single maximum at = for all 
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Figure 8: Contour plots of P (<pi;T) for (a) a/l = 1.0, (b) a/l = 1.3, (c) a// = 1.6 and (d) a/l = 1.9. 
Solid lines show contours at the labeled value. Dashed lines show maxima in P((pi;T). 




Figure 9: (a) Behavior of P(cj> = 0;T) with respect to a/l; to the left of the solid line, P{4> = 0) is 
a minimum, whilst to the right, P(cf> = 0) is a maximum. Directly on the line, P{<f> = 0) is a point 
of inflection. Dashed lines indicate the upper and lower bounds in a/l of the line, (b) The behavior 
of the maximum of P((f>i) T = 0) (that is, in the ground state) with respect to a/l. There is a kink 
in the line at a/l = 4/3, indicated by the dotted line. 



temperatures. The points separating the two regimes (max (P (<pi, T)) = and ^ 0) can be 
calculated analytically, as the points at which — — 



0. The calculated line in 



01=0 



overlap-temperature space is shown in Fig. [9] (a) . The upper limit of this line is the value of 
a/l at which (fi c (xi = x i+ i = h = a) = 0. For values of overlap equal to or larger than this, 
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Table II: Comparison of the densities of states for a square well tetramer chain, a /I = 1.0 and 
Act = 1.5 calculated in this work (LO^(Ek)/Al 4 , with the factor 1/4 as described in section ITTT1) 
and by Taylor (wf aylor \k), or in the original terminology) . Suppression of an unimportant 
multiplicative factor of 16ir 2 in the work of Taylor leads to the difference in the values; it can be 
seen that including this factor, the differences in the values are negligible. 



k u 5 {k)/4l A 


uf aylor \k) W4 (fc)/i 




U Iz.ybooyo 


O.OozUb 


1 AAAQQ/? 

i.OOOooo 


1 21.368822 


0.13531 


1.000071 


2 14.300928 


0.09057 


0.999908 


3 7.300283 


0.04626 


0.999342 


4 2.284838 


0.01447 


0.999924 


5 0.0633290 


0.004012 


0.999590 


6 0.035383 


0.0002222 


1.008395 



it is not possible for the polymer to exhibit 1-4 overlaps, and there is no steric hindrance to 
0i = states, which are the points of closest 1-4 approach. The lower limit of this region 
occurs at the point where the maximum of P{4>x] T = 0) (see Fig. [9] (b)) becomes zero. The 
value of a /I at this limit does not admit a simple interpretation or expression. 

B. Pentamer 

A comparison between the densities of states calculated here for pentamer tangent square 
well chains with those presented by Taylor is provided in Table [III Results are equivalent to 
four significant figures. The calculated densities of states as a function of a /I are shown in 
Fig. QUI Once again, we see the highest energy densities of states going to zero at a /I = 4/3 
as 1-3 interactions become "always on", combined with a kink in the density of states for 
the highest remaining energy. All densities of states tend to zero as a /I — > 2, where the 
available conformational space becomes zero. There are two further behaviors, not seen in 
the tetramer. The most obvious is that the density of the lowest energy state u 5 (6) becomes 
zero at a /I = a/7/2. For values of overlap larger than this, the pentamer has become so 
stiff that it cannot bend back on itself far enough to make 1-5 contacts. 
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Figure 10: Densities of states for pentamers plotted against a /I. 
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Figure 11: Concavity in the entropy and bimodal energy probability distribution for a /I = 1.55. 
(a) The density of states w^{k). (b) The entropy as a function of k; note the concavity at k = 5. (c) 
The probability distribution function P{k, T) at the "state coexistence" temperature. The function 
is bimodal, and the total weights of the two "states" (k < 5 (diagonal shading) and k > 5 (horizontal 
shading)) are equal. Dashed lines serve as a guide to the eye. 

A further interesting behavior is observed at intermediate values of a /I where the ground 
state oj§{Q) becomes the same order of magnitude as co>5(5). Indeed, for 1.53 < a /I < 1.56, 
1^5(6) > co> 5 (5). This gives rise to a concavity in the entropy S(k) = ks lnu(k) of the system 
with respect to energy at E = — 5e, which can be studied using the discrete analog to the 
second derivative, S"(k) = (S(k + 1) — 2S(k) + S(k — 1)); the function is concave if S"(k) is 
negative. The concavity results in a bimodal probability distribution function P(E;T,a/l) 
(illustrated in Fig. ( ITTT) ) . In analogy to the study of phase transitions, we find the line of 
temperatures at which the two peaks of these bimodal probability distributions have equal 
weight - a line of "state coexistence". This line is plotted alongside the data in Figs. [12] 
(energy) and [13] (heat capacity), and runs from a /I ~ 1.14 to a /I ~ 1.72. These end points 
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Figure 12: Ensemble average energies (E) /e for pentamers plotted against c/Z and temperature 
T. Solid lines show energy contours at the labeled value. The dashed line shows the "state co- 
existence" line. Note the discontinuities in the slope of the energy contours at a/l = 4/3 and 
cr/Z = a/7/2 (shown by dotted lines). 
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Figure 13: Configurational heat capacity Cy for pentamers plotted against a/l and temperature. 
Solid lines show contours at the labeled value. The dashed line shows the line of maxima in Cy with 
respect to temperature T. The dot-dashed line shows the "state coexistence" line — the position 
of the end points of this line near a contour line is purely coincidental. Note the discontinuities in 
the slope of the heat capacity contours at a/l = 4/3 and a/l = a/7/2 (shown by dotted lines). 



are at non-zero temperature, and occur where the curvature of the free energy at E — — 5e 
becomes zero. The end points are not associated with heat capacity divergences. 

The thermodynamic data shows the expected discontinuities in the slope of energy and 



heat capacity contour at the "magic numbers" a/l = 4/3 and a/l = y 7/2. The lower magic 
number corresponds to the loss of high energy states, as for the tetramer. The larger magic 
number, corresponding to the loss of the k — 6 state, gives a discontinuity in the energy at 
zero temperature (from E = — 6e to E — —he). As for the tetramer, bond length fluctuations 
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Figure 14: Ground state dihedral densities of states g§ (4>i, fa, 6) for (a) a/l = 1.0, (b) a/l = 1.3, 
and (c) a/l = 1.6, and (d) 55 (0i,02 5 5) for a/l = 1.9; the ground state is — 5e for a/l = 1.9. Lines 
show contours as denoted in figure legends. 

in real systems would act to smooth out these discontinuities in real systems. 

The pentamer also shows a line of heat capacity maxima, which lies at lower temperature 
than the "state coexistence" line. Both these lines show a discontinuity in slope at a/l = 4/3. 
Both lines are doubly reentrant, showing one maximum below a/l = 4/3, and another above 
a/l = 4/3. The line of heat capacity maxima connects with the discontinuity in energy at 
a/l = v/7/2. 

The dihedral behavior of the pentamer at zero temperature (ground state) is shown 
in Fig. [HI in four representative plots of the ground state dihedral densities of states 
95 (01? 02; 6) and c?5 02? 5) . These are equivalent to unnormalized dihedral probability 
distributions for the system at T = 0. We see that, for the three lowest values of a/l, the 
probabilities are peaked at points on the 0i = 02 diagonal; the dihedrals prefer to take the 



same sign. This continues to the "magic number" a/l = y (3 + v / 5)/2, where, as for the 
tetramer, 1-4 overlaps can no longer occur, and steric effects no longer prevent cis conforma- 
tions. As the degree of overlap tends to this number, and the amount of steric interference 
decreases, the maxima move closer to 0i = 2 = 0. For a/l > J (3 + \/5)/2, the probability 
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distributions become unimodal at 4>i — 4>2 = 0. This should be compared with the behavior 
for the tetramer (see Fig. [9]), where (fi ma x(T = 0) is zero for a /I > 1.48; this effect is due to 
the additional steric interference from 1-5 overlaps. 



V. DISCUSSION AND CONCLUSIONS 



In the previous section, it has been shown that the tetramer and pentamer show a rich 
and surprising range of behaviors. Specifically, these are "magic numbers" of the overlap a/l 
where the derivatives of the densities of states change discontinuously, maxima in specific 
heat with respect to temperature, and a region of bimodal energy probability distributions, 
reminiscent of a first-order transition in bulk systems. In general, the behavior of long 
polymer chains cannot be directly inferred from the behavior of very short chains such as 
those studied in this work. If, however, interactions between monomers widely spaced along 
a chain can be neglected, the behavior of very short chains can be used as a basis for a spin 
chain model. Such interactions may be neglected when chains become very stiff (at, for e.g., 
large values of a //, or after helix formation). In this case, the behavior of the very short 
chains may be considered the "building block" for the behavior of longer chains. 

The "magic numbers" which are observed correspond to discontinuous changes in the 
derivatives of the densities of states. At a/l = 4/3, 1-3 contacts become "always on" and 



high energy densities of states become zero. At a/l = y (3 + \/5)/2, the chain becomes so 
stiff that it cannot bend back upon itself far enough for 1-4 overlaps to occur. Similarly, 



at a/l = a/7/2, the chain becomes so stiff that 1-5 contacts can no longer occur, and 
the ground state for the pentamer is lost. These discontinuities in the densities of states 
are associated with discontinuities in the energy and compressibility with respect to the 
parameter a. The "magic numbers" are similar in principle to the "cut-off" A values noted 



by Taylor 



12j for tangent chains — this work has not examined the effects of changing 



the well width parameter A, but it is obvious that the values of these "magic numbers" 
will depend upon that parameter, and that "cut-off" values of A will also exist for this 
model. As has been noted above, the discontinuities across lines of constant a/l in this 
system will be smoothed in simulations with variable bond length, however, the effects 
should still be visible. We particularly note the sudden loss of stability of the "helix 1" phase 
at a/l ~ 1.675 in previous simulation work {3]] (see Fig. [T]). Given the 10% bond length 
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fluctuation allowed in those simulations, this loss of stability may coincide with the magic 



number at a/l = y (3 + y/E)/2 ~ 1.618, suggesting that the more tightly wound "helix 1" 
phase is stabilized by steric interference of 1-4 contacts. This supposition is supported by the 
observed loss of double-peaked dihedral angle probability distributions for overlaps above 
this "magic number", suggesting that the more loosely wound "helix 2" phase is stabilized 
by steric interference between monomers spaced further along the chain. 

The low temperature maxima in the specific heat for these short polymers appear to be 
a continuation of the specific heat maxima observed for short tangent chains Q|- We 
follow these previous works in interpreting these maxima as signatures of collapse to close- 
packed, low energy conformations. This interpretation appears confirmed by the presence of 
bimodal energy probability distributions for the pentamer, with a line of "state coexistence" 
which roughly parallels the line of maxima. 

For the tetramer, the line of specific heat maxima shows re-entrance below a/l = 4/3, 
having a maximum with respect to temperature. For the pentamer, both the line of specific 
heat maxima and of "state coexistence" are doubly reentrant, showing maxima below and 
above a /I = 4/3. The re-entrance of the "state coexistence" line can be easily explained by 
reference to the densities of states shown in Fig. [101 Consider the system for a/l > 4/3. 
For overlaps just above this point, the ground state density of states (the entropy of the low 
energy state) is increasing while all other densities of states are decreasing with increasing 
overlap. Hence, the low energy state becomes more stable, and coexistence moves to higher 
temperature. The ground state density of states soon begins to decrease, but as long as it 
is decreasing more slowly than the higher energy density of states, its stability continues to 
increase. However, on closer approach to a/l = \/7/2, the ground state density of states 
decreases faster than the higher energy densities of states, and stability decreases. The 
same argument holds for the line when a/l < 4/3. If we interpret the maximum in heat 
capacity as a result of structural competition between the ground state and higher energy 
states (following Stanley, et al. [19|]), we can make the same argument for the re-entrance in 
the lines of maxima for both the tetramer and pentamer. Physically, increasing the overlap 
of the chain makes configurations with lower energy (more contacts) more likely at first (as 
monomers are "drawn into" each other's square wells), but then begins to cut into these low 
energy states as the chain becomes too stiff to bend back upon itself and make contacts. 
We attribute the re-entrance of the stability of the "helix 1" phase in previous work to this 
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same competition between effects. 

Though it seems reasonable to attribute the behavior of the phase boundary between 
the "helix 1" and globule phases to effects seen in the pentamer, it should be noted that 
the state coexistence seen in the pentamer is represents collapse of the pentamer, rather 
than helix formation. Though the dihedral probability distributions shown in Fig. [Q] do 
show double peaks at non-zero dihedral angles, this is not a sufficient criterion for helicity. 
The cross-correlation coefficient of these distributions is not significantly above zero; the 
total statistical weight associated with dihedrals away from the peaks is still large enough 
to outweigh the correlated peaks. However, the clear double peaked structure does suggest 
that the physics necessary for helix formation is contained in these simple, small systems, 
particularly in the steric interference due to 1-4 overlaps. 

While these results appear to clarify certain behaviors observed in simulations, they do 
raise further questions. Under the interpretation we have offered here, the nature of the 
"helix-2" phase is unclear; this phase is observed to be stable up to all = 1.9 in simulation 
[3J], where the chain is too stiff for 1-5 overlaps to be the root of the observed chirality. 
Further, the question of how the helix transition connects (or does not connect) to the 
crystallization-like transition observed in simulations for the tangent chain system remains 
unresolved. Follow-up work, developing a spin chain model for helix formation using the 
results presented here, is underway; it is hoped that this approach will shed light upon these 
questions. 
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